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A new dynamical cascade code for decaying hot nuclei is proposed and specially adapted to the 
synthesis of super-heavy nuclei. For such a case, the interesting channel is the tiny fraction that will 
decay through particles emission, thus the code avoids classical Monte-Carlo methods and proposes 
a new numerical scheme. The time dependence is explicitely taken into account in order to cope with 
the fact that fission decay rate might not be constant. The code allows to evaluate both statistical 
and dynamical observables. Results are successfully compared to experimental data. 
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I. INTRODUCTION 



Since the pioneering work of Bateman [1] who developed a general equation for radioactive decay chains, the model 
has found many applications in nuclear physics and in other fields. Purely statistical version of the problem, when 
the time is integrated out, has been extended to multi-channel cascade process of decay of exited compound nuclei 
that can evaporate neutrons, protons, alphas. . . and successfully compared to experimental data. The original time 
dependent Bateman equations can also easily be extended to the multi-channel process, provided that the decay widths 
are time- independent. However, the decay width of the fission channel is now known not to be time-independent, 
so Bateman problem should be solved numerically. Actually in fission of excited compound nuclei, it is known that 
it takes a finite time to build up the quasi-stationary probability flow over the fission barrier. Therefore, the time 
evolution should be taken into account. For recent reviews on the dynamics of fission, see Refs. [2,3]. 

One of our purposes is to apply the method to evaluate the production cross section of the superheavy elements. In 
such reactions, most of the compound nuclei will undergo fission, and thus we need to evaluate the very tiny fraction 
that will decay only through particles emission. Therefore, classical Monte-Carlo methods [4,5] are not suitable. We 
propose here another numerical scheme. Furthermore we construct a new statistical code named KEWPIE (Kyoto 
Evaporation Width calculation Program with time Evolution) . 

Fusion reactions that lead to the production of super-heavy elements involve often a combination of rather heavy 
ions for the entrance channel where it is well known that fusion is strongly hindered. The extra push parametrisation 
gives only a qualitative behaviour of this fusion hindrance, so more precise calculation should be done. There are 
several way of taking into account the hindrance, one of which is to use the two step model [6-8]. The code is designed 
to accommodate fusion probabilities calculated by such a model, in addition to the probabilities given as usual by 
transmission coefficients whose calculation with the full proximity potential is included in the code as standard. 

II. DYNAMICAL FRAMEWORK 



A. Single chain 



1. Populations as a function of time 

Before generalising the method to the multi-channel decay problem, we first recall Bateman's results considering a 
single chain starting from a single nucleus. The equations then read, 

dp 

—l^^r\Pi with Pi(o) = i, (1) 

dp 

^TiPi-rlP2 with P2(0) = (2) 

dP 

—1^T2P2'T\P^ with P3(0)-0 (3) 
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dP 

^ = r„_iP„_i - r*„p„ with p„(o) = o, (4) 

where is the total decay width for the nucleus n and r„ is the evaporation width, whatever the particle is. Pn{t) 
denotes the population of the nucleus n at time t. If the decaying widths are time independent, it is very easy to solve 
analytically these coupled differential equations, by taking the Laplace transform as Bateman did. This leads to 

PM - (5) 



where P(s) is the Laplace transform of P{t). The solution is then obtained by taking the inverse Laplace transform, 

Pi{t) = e-^'^' (8) 
P2W = ^^(e-^^*-e-r^*) (9) 

^ 2 ~ ^ 1 



fe=l i=l IWiV^J 



(10) 



where we have supposed that the F^'s are all different from each other. A special attention should be drawn on the 
last nucleus of the chain, for which the solution might differ from the above when F* = 0. For such a case, 

'•'max 

r+oo 

f"n^a. = r„_.-l / dtPn^^^.,{t) (11) 

Jo 

= n E rm =. (r*-F*)- ^^^^ 

fe=l i=l »llj#»V 3 il 

From the populations, we can then calculate any observable. This solution is only valid for time independent F*'s 
and then for a given excitation energy for each nucleus. It is then not always applicable, but it will be very useful to 
test our numerical scheme. 

2. Statistical observables 

Assuming that only neutrons are evaporated along the chain, the total width for the isotope n that can also undergo 
fission, is then F^ = F^^ + F„, where F^^ corresponds to the fission width. The probability to emit exactly n neutrons 
prior fission is 



Pn-l 



dtriP^_,{t) (13) 







=&n& (14) 



F* J--'- F* 

" i=i » 



which is in agreement with the reults of Ref. [9] . For the end of the chain. 



Pn^..= ""f[ (15) 
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These results are commonly used in statistical decay codes. The average pre-scission neutron multiplicity is simply, 

"ma x 

(n)= ^(n-l)p„. (16) 

n=l 

The total number of nuclei that fission is that TV/ = 1 if r„^^^ ^ and Nf = 1 — Pnmax ^l^^- 

3. Fission time 

The average fission time is an interesting observable to test the dynamics of the cascade because it can be measured 
directly in some cases. We define it as 



+°° dPt 



t^^dt (17) 



-I 

^ / tTfMt)dt (18) 



B. Multichannel scheme 

Excited compound nuclei can evaporate other particles than neutrons. Protons, alphas should also be taken into 

consideration. The populations will be labelled as Pi_j with i being the number of evaporated neutrons and j the 
number of evaporated protons. In order to keep a triangular shape of the matrix, the populations will be ordered 
following the number of evaporated nucleons, starting with the neutrons. Or, in a more formal way, n = + j 

where n is a ranking number for each nucleus produced in the decay process. 

Then, similar differential equations can be written, taking into account, neutrons, protons, alpha evaporation and 
fission. The first ones read 

(19) 
(20) 
(21) 
(22) 
(23) 
(24) 



dPo,o 
dt 


= — r^Po,o 


with Po,o(0) = 


1, 




dPifl 
dt 


= oPo,o 


~ o-Pi,o 


with 


Pi AO) = 





dPo,i 
dt 


= ^0,0-^0,0 


~ r'o.i-Po,! 


with 


Po,m = 





dP2fl 

dt 




~ r2^QP2,0 


with 


^2,0 (0) = 





dPi,i 
dt 


= r'o,i-Po,i 


+ r?,oA,o- 


- r*i_iPi,i with 


Po,i(0) = 


dPo,2 
dt 


= ^0,1-^0,1 




with 


^^0,2(0) = 






^ =r^,1^2,l+r?,2Pi,2 + r^,0^0,0 -4,2^2,2 with P2,2(0)=0 (25) 



Here r*^- corresponds to the total decay width of the nucleus ij, including fission. The other F's correspond to 
evaporation widths. The Laplace transform of the first equations is very similar to what was done with a single chain, 
up to Eq. (23) for which some differences appear: 

pP pP pn 
p / \ _ 0.1 0.0 1,0 0,0 

" (ri.i + .)(4,i + s){rl, + s)^ [Ti, + s){T% + ^)(rU + ^) ' ^ ' 

The result corresponds to the sum of single chains over all possible paths. Just considering neutrons and protons 
evaporation, there are ^^j^ possible paths from the nucleus (0,0) to the nucleus Then, taking the inverse 
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Laplace transform, one ends up with a general formula that is of course the sum of single chain terms over all possible 
paths. This general property due to the linearity of the equations is also true for the observables, for which time 
integration should be carried out. It is possible to recollect similar terms together, but it will make formulae more 
complicated. 

The complexity of the problem is eventually due to the multi-channel scheme as in any cascade code, not to 
the dynamics that could be easily added exactly to a statistical code, assuming that the F's are time independent. 
Nevertheless, these equations have the merit to be exact, but based on an approximate average model. Since we are 
not only interested in the main stream decay channel but also small decay paths, a numerical solution is necessary. In 
addition, this model docs take into consideration 7 emissions that cool down nuclei. These results will then be useful 
for testing numerical scheme and understand physical results. 



C. Numerical scheme 



For the sake of simplicity, we will only describe the scheme for a single chain. It could be easily extended to a 
multi-channel cascade. In the code each time step is larger than the previous one: (5ij_|_i = a6ti. The alpha parameter 
permits to have an increasing time step along the cascade, and thus permits to take into account processes that can 
have very different time scales. The parameter a is usually set to 1.2 which permits to cover a wide range of excitation 
energy. Larger values of the parameter a may induce a divergence in the fission time calculation (equation 18) due 
to the multiplication of a small probability by a long time. The initial value Sto is set to ITi/MeV. This choice is 
made so that initial time scale is longer than the mean collision time (around O.lTi/MeV) and shorter than the mean 
evaporation time (few tens of h/MeV). This value is physically reasonable for a wide range of application. 

During a time Sti, the population of a given nucleus Pn+i increases by the decay of a mother nucleus (if exits). At 
a given excitation energy i?*, this feeding term is Pn{ti)N^ with 

7V^ = ^^(l-e-^"**0- (27) 

During the same time this population is decreasing according to an exponential law : P„_|_i(tj -1- Sti) = Pn+i{ti)Dn-\.i 
with -D^_|_i — e~'""+i***. Finally the iterative equation for the population P„ reads, 

Pn+l{ti + Sti) = Pn+l{ti)Dl+^ + Pn{ti)Nl (28) 

This is for an averaged excitation energy. 

In the multichannel scheme the representation is more complicated due to the interconnection between the different 

ways of disintegrations. However, the problem can be treated in the same framework. We can calculate the number 

of nuclei going to fission at each step by looking at the quantity Pn-i{ti)F^^ with F^_i = ^(1 — e^'""''*') for each 
nucleus n at the time ti. The fission time can easily be deduced with a discretised version of equation (18). With 
this framework we can solve more sophisticated dynamics including cooling by gamma emission or a transient time 
during which fission probability is reduced. 



D. Energy spectra 

One of the particularity of the KEWPIE code is that the energy spectra of the producxid nuclei arc completely 
calculated and processed. Monte-Carlo methods were not considered because we are interested in the tiny fraction of 
events leading to super- heavy elements. For a given mother nucleus n with a given excitation energy E^, the energy 
spectrum of the daughter nucleus n -|- 1 is proportionnal to jiE* — > E*_^-^) entering the evaporation width. 



n{E*J= dej{E*^^E*^+,), 

JVcoul 



(29) 



where e = i?* — -Bj^+i is the energy taken by the evaporated particle and Vcoui is the Coulomb barrier between the 
daughter nucleus 4and the evaporated particle. The former is naturally zero for neutrons. 7(-E'^ ^n+i) is the 
partial width for a given value of e, 

7(K-K+i)cx^^^ (30) 
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where p represent the level of density. The detailed formula is given in section III A. When the mother nucleus has 
itself an energy spectrum Sn{E*,t) normalized to the population, J^°° Sn{E*,t) dE^ = Pn(t), the energy spectrum 
of the daughter is feeded by 

SSn+i{E:^„U + 6U)= / dE*^Sn{E:,u y'- " "+i^ (i-e-^^^*-)- (31) 

Jo J- n 

In this equation, the total width is already integrated over the excitation energy i?*. The spectrum of the daughter 
nucleus will be used as an input for the next step of the cascade. The spectrum of the mother nucleus is modified with 
respect to the decay toward the daugther. Practically, in the code, the energy spectrum is evaluated by discretizing 
the previous integral in A'' excitation energy bins. 



III. PHYSICAL INGREDIENTS 



The physics of the problem is then in the various widths entering in the previous equations. In this part we will 
specify our choice for super-heavy elements. In the code the integrals are evaluated by discretizing the sum. 



A. Particle evaporation width 

Evaporation process in excited compound nuclei was modelised a very long time ago by Weisskopft [10] with the 
principle of the detailed balance. But this modelisation does not take into accoimt the angular momentum explicitly 
though it is a very sensitive parameter to modelise the competition between fusion and evaporation. We prefer the 
Hauser-Feschbach approach [11] which gives the width of the disintegration of a compound nucleus C at excitation 
energy E^ and spin J toward a nucleus B by the emission of a particle b of spin 7^, 



1 pE^-Bl, J" b+-f -B J+S b 



(32) 



In this equation represents the level density at the ground state and ^^5^^ the transmission coefficient. Detailed 
calculations of those coefficients are described in section IIIG and IIIJ, respectively. The binding energy of the 
particle h in the nucleus B is represented by B5. The variable Sh is the kinetic energy of the evaporated particle b. To 
obtain the total width, a summation over all the kinetic energies of the emitted particles is done from the Coulomb 
barrier between the nucleus B and 6, named VcouU to the maximal value reachable E* — B},. 

We will use in the remaining part of this publication a simpler notation F^ in place of ^'^c-j-e' }^{b-ii,} where ^ is 
the evaporated particle. 



B. Fission width 



In order to calculate the fission width, we need to know the height of the fission barrier Bj which is calculated in 
detailed in section IIIF. Then the width is given by the Bohr formula [12], in which a transmission coefficient factor 
Tf{ef) is included. 



1 r'^c-"s 

rf(^7;B£}^{/;J;B,} = 2wp3{E^;J) J_b^ defp\E*c - Bf - Sf, J)Tf{ef). (33) 

In this equation and represent the level densities of the nucleus at the ground and saddle points, respectively, 
see section IIIG. The transmission coefficient is given by the Hill and Wheeler formula [13]. This factor takes into 
account the effect of sub-barrier fission, 

Tfisf) = (34) 

1-hel ^) 
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where the parameter LUsd is related to the curvature of the fission barrier at the saddle point. We can define in the 
same way LUgs as the curvature of this potential for the deformation around the ground state. In this code hugs and 
hiOsd are usually set to IMeV. The variables ojsd and cOgs are also useful in order to add the Kramers [14] and the 
Strutinsky [15] correction factors to the fission width : 

r{C;J;E*}^{/;J;E/} = ^^^{C^J;E'^}^{f;J;Ef}- (35) 

The Strutinsky factor <S is given by 

S=^. (36) 

This pre-exponential factor includes temperature T {E* = aT^) which cancels out the temperature dependence 
implicitly included in Bohr- Wheeler formula. 
The Kramers factor IC reads, 



K = yjl + xl - Xk, (37) 

with 

= (38) 

Those formulae take into account the effect of the viscosity of the nuclear matter against the fission process, where 
/?/ is the reduced friction coefficient i.e. the friction coefficient divided by the inertia mass. Its value is not well known 
yet and set around lO^^s"^ to lO^^s"^. 

The notation for the width given in equation (35) will be simplified as Tj. 



C. Gamma emission width 



Another mode of de-excitation of the compound nucleus is the gamma-ray emission. This process reduces the exci- 
tation energy of the nucleus without changing the nucleus itself. That de-excitation can lead the nucleus energetically 
under the other thresholds of disintegrations (fission and evaporation) and so contributes to the stabilisation of the 
excited compound nucleus. The width of the gamma emission is governed by the following formula [16,17] 

^{C;J;E',}^j= 2^o(E*-J) ^ ^ / ds,pc{E*c - e,; 1)^6^+^ (39) 

In the calculation code KEWPIE we only consider the transitions with L=l and L=2. Integration is done between 
a minimal -E^^ and maximal E^j^ allowed for the energy of gammas : E}^^ = 1 MeV and E}^j^ = 7 MeV for L = 1 
and E^^ = 1 MeV and E^j^^ = 4 MeV for L = 2. The factors and ^2 are parametrised as follow 

^1 = 4.24 lO^^^t (40) 
C2 = 3.00 10-^2^t. (41) 

Theses values only depend on the mass of the nucleus A. 



D. Mass and Shell correction 



Values of nuclear mass Ma and shell corrections AEsheii of nuclei are taken from the M0ller et al mass table [18]. 
In this table are referenced both experimental and theoretical values for mass and shell correction. If the experimental 
value is known, this value is taken. In the other cases the masses calculated by the M0ller are used, which include 
macroscopic liquid drop energies and the shell correction energies obtained by Strutinsky method. Another alternative 
is to use the mass table proposed by Koura [19]. 



6 



E. Fissility 



The fissility of nucleus is used as an index to evaluate many quantities. It is defined by the Coulomb energy divided 

by twice of the surface energy of the nucleus and can be given practically in several ways. The simplest standard 
definition is a; = Z'^/5QA, with A and Z representing the mass and charge numbers of the nucleus, respectively. For 
heavy nuclei, we will prefer another parametrisation taking into account the effect of isospin. The ratio between 
neutrons and protons changes the sensitivity of the nuclei to the fission process. More over this formula of the fissility 
is deeply related to an evaluation of the fission barrier that is described in the next section. A formula well adapted 
to massive systems is given by K. H. Schmidt et al. [20] which is used in the present code, 

X = - (AO) 

49.22A(1 - 0.3803/2) _ 20.4897^) ' ^ ^ 

where I = N — Z = A — 2Z is the isospin. We principally use this formula that, as mentioned before, is used for the 
parametrisation of the liquide drop part of the fission barrier III F in order to keep a coherence between the various 
parts of the code. 



F. Fission barrier 



There are many ways to parametrise the fission barrier Bfo, the Liquid Drop Model (LDM) fission barrier for J = 0. 
The code is basically designed to use two parametrisations. The first one is the classical Cohen-Plassil-Swiatecki (CPS) 
[21] parametrisation that is multiparametric fit of the height of the fission barrier as a function of deformation of the 
nucleus based on a liquid drop model. 

As the purpose of the KEWPIE code is mainly to deal with heavy and super-heavy elements, another parametri- 
sation is also implemented in the code. In the paper related to heavy elements, K.H. Schmidt al [20] proposed a 
phenomenological prametrisation of the liquid drop model fission barrier as a function of the fissility. This parametri- 
sation is obtained by adjusting parameters with data on a wide range of heavy nuclei compiled by Dahlinger et al. [22]. 
The results of the fitting permit to reproduce the data related to the disintegration of the heavy elements relatively 
well over a wide range of nuclei. Then with this parametrisation, in the code, the fission barrier is taken as, 

_ 0.7332 10-^2 

^fo = —1 , (43) 

A3 

with u = 0.368 — 5.057x + S.dSx"^ — 8.71a;'^ where x is the fissility of the nucleus under consideration. The formula 
represents the height of the barrier for a null angular momentum and without taking into account the shell correction 

energy AE shell- Thus it has to be modified to take in account those effects. First, the evolution from the spherical 
shape at the ground state to the deformed state at the saddle point leads to a modification of the barrier due to the 
change of rotational energy: 

Bf, = Bfo + {ETor^ - EZr"), (44) 

where E^^^'^'''^ and E^^^^"^ represent the rotational energies at the saddle point and the ground state, respectively. 

Subsequently the shell energy correction AE shell mentioned in section IIIF is substacted to obtain the true value 
of the height of the fission barrier, as 

Bf=Bf^-AEsheU- (45) 

We are assuming that the shell correction at the saddle point can be neglected. This is the final value of the fission 
barrier Bf that is used in the calculation of the fission width HID. 



G. Level density 

In the program KEWPIE the evaporation of particles is calculated for each partial cross section as a function of 
the angular momentum J. So we need a level density function with angular momentum specified and we take one 
given in the textbook by Bohr and Mottelson [23], 
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^1 , ^ . 2 2./a(E'' -Erot) 

piE\J) = -6-H2J+l)a(-) (46) 
48 \2'-^J {a{E* - Erot))^ 

The rotational energy of the nucleus is given by 

_ J{J+l)h' 

™* - 23 ■ *- 

In this formula J represents the angular momentum quantum number of the nucleus. To obtain the rotational energy 
we need to calculate the moment of inertia 9 which is detailed in the next section. 

The level density is the very key physical quantity in the calculations because it comes in all the disintegration 
widths, whatever the processes are, i.e. , evaporation, fission and gamma emission. 

H. Moment of inertia 

In the program KEWPIE we will assume that the nucleus at the ground state is spherical. It is reasonable for 
superheavy elements in the sense that the most stable nuclei that we hope to form are magical nuclei and neighbours. 
Since their fissility is high, the nuclei are not very deformed also at the saddle, so the evolution from the ground state 
to the saddle point is rather small. The moment of inertia for the spherical nucleus is evaluated as a rigid body [24] 
as follow, 



Aground • (^^) 

on 



The moment of inertia at the saddle point takes into account the deformation of the nucleus by multiplying the 
momentum at the ground state by a modification factor given by Hasse and Myers [25] , 

^saddle = "Aground ( 1 + -^V ( 1 + -^^V] ) • (49) 



6" V 255 _ 

In this equation y = 1 — x with x the fissility of the nucleus. We have limited this expansion to the second order in 
the code, as it is suggested in reference [25]. Keeping this limitation in equation (49) at the second order in y we can 
express "^saddle as a second order function of a2 (factor of the second order Legendre polynomial development of the 
nucleus shape) by the following formula, 

cv ... -o .^1^1 ^424 , 



^saddle = Aground + + "^^2 j • (50) 

For the level of density at the ground state or the saddle point we will use the corresponding parametrisation for 
moment of inertia. The saddle point deformation can be determined by shell correction energy. In that case, it can 
be read from a file in order to accommodate results of other microscopic calculations. 

I. Level density parameter 

In this program the level density parameter a is given by the formula of Toke and Swiatecki [26], 

-(1 + 3.114^+5.626^), (51) 



14.6r Ai Ai 

where Fi and F2 arc related to the deformation of the nucleus as the variation of the surface energy and curvature 
energy, respectively. At the ground state, F2 and F3 are equal to 1, as we assume no deformation of the nucleus. To 
take into account the effects of deformation at the saddle point, we use again the parametrisation done by Hasse and 
Meyer [25], 

= ! + (52) 



45" V 765" 
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where we have again hmitcd the expansion to the second order. In those relations y = 1 — x which is expected to be 
very smaU for the heavy or super heavy nucleus. 

At the ground state, an effect of the energy-dependence of the shell correction is taken into account, according to 
Ignatyuk's prescription [27], 

aground = a{l + f{E*)^^^), (54) 

with 

f{E*)=l-e~''^. (55) 

In the previous equation E* is the excitation energy. The parameter Edamp has been adjusted with measured level 
densities and represent how much the shell correction is damped by the excitation energy. Its value in the code is 
taken as 18.5 MeV. At the saddle point we assume no shell correction. 



J. Transmission coefficient 



Assuming as usual that, 7^j^g^(6, £&) w 7;^ (&,£(,) = Ti, we have to solve the Schrodinger equation, 

1(1 + 1) 



P [Vcouiir) + VNucir)]^ V(r-) = 



(56) 



from r = to infinity in order to exploit the asymptotical properties of the wave function. At very large values of r 
we can calculate the Si matrix and then obtain the transmission coefficient Ti = 1 — \Si\'^ . In this equation Vcoul and 
Vnuc represent coulomb and nuclear potentials, respectively. 

At infinity for a given angular momentum I the general solution of the Schrodinger equation reads 

V'(r) (X n\~\r],kr) - Sn\'^\ri,kr), (57) 

with 



'l-6{^\ri, kr) = Qiiri, kr) + iTi{ri, kr), 
'H\~\ri, kr) = Qi{ri, kr) - iTi{ri, kr). 



(58) 
(59) 



In those equations !Fi{ri, kr) and Qi{ri, kr) represent the regular and irregular Coulomb wave functions, respectively. 
By using the values of the function i]j{r) at two closes values r\ and ro with ri > ro » we can obtain the value of 
the Si matrix by the following formula, 



5. = ^' 



(-) 



[t], kr o)ip [r i) -H\ \ri, kri)-tjj{ra) 



n\'^\riMo)^{ri) 



H[+^(??,fcri)V(ro)' 



Then we can easily obtain the value of the transmission coefficient Ti. 



(60) 





Neutrons 


Charged particles 


Ur 


47 - 0.27^teb 


54 - Q22Eiab + 0.4 A + 24^^^ 


Rr 


(1.332 - 7.610-*A + 410-^^2 - d,lQ-^A^)A^ 


1.17A3 


ttR 


0.66 


0.75 


Uso 


7 


6.2 


Rso 


(1.266 - 3.710-^^^ + 210-6^2 _ 410-9^3)^^ 


1.01^^ 


aso 


0.48 


0.75 


Wi 


0.22Eiab - 2.7 


Q.22Eiab - 2.7 


Ri 


(1.266 - 3.710-M + 210-*'^2 _ 41(3-9 A^)^* 


1.32-4 IT 


ai 


0.48 


O.T)! + 0.7'",-^ 


WsF 


9.52 - Qm'iEiab 


11.8-0.25S,afc + 12^i^ 


RsF 


(1.266 - 3.710-4^ + 21Q-^A^ - 410-^3)^* 


1.32A3 


asF 


0.48 


0.51 + 0.7^-/^ 



TABLE I. Values of the parameters entering the optical nuclear potential 
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To solve this problem in the KEWPIE code we make all the calculation numerically, using the most efficient 
numerical method for evaluation of the wave function, namely the modified method of Noumerov [28] . 

The nuclear optical potential depends on the type of the evaporated particle. Its real and imaginary parts are given 
as follows, respectively, 

rr ^^ tt ft n ^ , UsoJ>^l df{r, Rso,aso) 

UNucir) = URf[r, Rr, an) H , (61) 

r ar 

WNucir) = -Wif{r,Ri,ai) + Wsr^ai^^^^^^^j^^^, (62) 

with 

f{r,R,a) = (63) 

1 + e^ 

In these equations J represents the angular momentum and A,r the pion Compton wavelength. The other parameters 
are listed in table I. We use the values given in Ref. [29] and Ref. [30] for neutrons and charged particles, respectively. 

Those parameters have been obtained by fitting the angular distributions of proton and neutron scattering. There 
are no common parametrisation for both neutron and protons. We assume that all charged particles have the same 
behaviour as protons, with respect to optical potential. 



IV. VALIDATION AND APPLICATIONS 



A. Statistical observables 



In order to test the code, its results are compared to data coming from various experiments related to the production 

of heavy and super-heavy elements. Those comparisons are done in two steps. First we will look at the results of the 
evaporation process from a fusion cross section calculated with a full proximity potential, where the fusion is supposed 
not to be affected by fusion hindrance. Next, we will look at the evaporation residue cross sections of some superheavy 
elements, where the fusion is calculated by the two step model which takes into account the fusion hindrance in a 
realistic way. 

We look at the xn residue cross section (residue after evaporation of x neutrons) for various combinations of 
projectiles and targets chosen so as to form the same compound nucleus ilgiVo^"^. In fig. 1 are represented the 
cross sections of evaporation residues of the xn reactions as a function of the excitation energy for the ^^C +^^^ Cm, 
16q _|_242 22jYg _|_236 jj ^^^^ "^^ M g +^''^ Th reactions. The experimental data are taken from the reference [31 33] 
for the +246 Cm , +^42 Pu and "^"^Ne +^36 U and "^^Mg +^32 Th reactions, respectively. In this figure the 
symbols represent the experimental data and the lines the results obtained by the code KEWPIE. 

In the present case, as stated above, the fusion probability and thus fusion cross sections, are calculated by using a 
full proximity potential [34] as nuclear attractive potential. The fusion hindrance is not taken into account, because 
the products of the charges of the projectiles and targets are lower than 1500 — 2000, which is the empirical value 
above which the fusion hindrance appears [35,36]. The calculations of the fusion probability are the same as the 
corresponding part of the HIVAP code [37] . 

KEWPIE code is designed in such a way that there are as few free parameters as possible. The only unconstrained 
value is the friction coefficient (3 for the Kramers factor (see equation (37)). The value of (3 is taken as 5 lO^^s"^. 
Setting the /3 to 5 10^^ (consistent with the one body wall and window formula) reduces too much the fission and 
this lead to an increase of the residue cross section by roughly one order of magnitude which obviously overestimates 
the data as shown in figure 2. No other free parameters are used. It might suggest that the fusion hindrance exists 
somehow in those systems. 

The code takes into account the competition between the evaporation of neutrons, protons and alpha particles 
as well as gamma emission, and thus a large number of paths in the disintegration process. However, the charged 
particle evaporations can be neglected without changing drastically the results. The physical justification comes from 
the huge Coulomb barrier that makes the probability of charged particles emissions very low. That approximation 
could be used to save calculation time. 
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FIG. 1. Cross sections of the residues coming from the 4n,5n and 6n disintegration processes for C+ Cm, 0+ Pu, 
'^'^Ne+^^^ U and ^^Mg+^^^ Th reactions. The symbols are representing the experimental data and the curves are the results of 
the calculations with the KEWPIE code. 

The figure 1 shows that an overall global agreement is obtained even with the strict conditions we have put and 
without adjusting free parameters that arc usually employed in statistical analyses of data. 

To complete the evaluation of the code we looked at a reaction with a more symmetric entrance channel to produce 
nobelium isotopes: ^*Co +^°^ Pb ^154 No^°'^. The product of the charges of the target and the projectile {Zp x Zt = 
1640) indicates that the fusion process described by the transmission of the barrier with the full proximity potential 
is still enough to give a good evaluation of the fusion cross-section. The residue cross section of the xn evaporation 
processes are shown on figure 2. 

As in the previous calculation shown in figure 1, /3 is set to 5.10^°s~^. As we can see in the figure 2 the calculation 
well reproduces the experimental data as a whole. If we look more in details, we can see that the peaks of the 
distribution are always located at the suitable positions. The main discrepancies come from the widths of the peaks. 
The calculated distributions are a little narrower than those of experiment. This small disagreements would come 
both from the experimental and the theoretical accuracies. In experiment there is a kind of averaging over a range 
of excitation energy of the residue cross sections due to the resolution of the detection set-up and the lose of energy 
of the beam in the target, those effects making the distributions broader. On the theoretical side, the problem of 
the width of the evaporation can come from the assumption made on the duration of one decay step. In the usual 
calculation code the time between two steps is suppose to be infinite. But if we caclulate each step with a finite time, 
the distributions become wider due to the progressive feeding of some nuclei. Moreover the inclusion of a transient 
time should especially affect the high energy part of the residue spectrum and make the whole distribution wider. 
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FIG. 2. Cross sections for the fusion and In, 2n, 3n and 4n evaporation processes from the reaction *^Ca +^°^ Ph. The 
symbols represent the experimental data and the curves the results of the KEWPIE code. 

We can conclude that in the reaction of production of the nobehum isotopes the fusion cross sections calculated 
with the full proximity potential give rise to an overall good agreement with the data, combined with the present 

statistical decay calculation. 

B. Dynamical observable 

In order to test the dynamical description of the code, a suitable observable would be the fission time which could 
be easily calculated by using equation (18). The case of ^ssjj jg particularly interesting because its fission barrier has 
roughly the same value as the neutron binding energy. Therefore, competition between the two channels occurs all 
along the disintegration chain and gives rise to very long fission times. The fission times were measured directly by 
crystal blocking techniques which are nuclear model independent [39] . 

Figure 3 shows a comparison between theoretical calculation done with KEWPIE code and data from reference 
[39]. We take only into account the competition between fission and neuton evaporation. Moreover, the emission of 
gamma ray are taken into account in this process since its cooling effect at the end of the decay, when low exitation 
energies are reached, affects the fission time. Those approximations are done in order to save calculation time but 
are resonable with respect to the evaluation of the fission time. In this experiment the angular momentum of the 
fissioning 238 [/ jg suppose to be lower than 10?l. In the calculation three values of the angular momentum Ih, 5h and 
9h are taken. 

The upper part of the figure 3 shows the fission time (define by equation 18) as a function of the excitation energy. 
By a careful study of the output of the calculation, we noticed that the average value of the fission time is affected 
by very slow fission events at the end of the decay chain. The calculated results show some oscillations that are 
damped at high energy. The origin of these oscillations is mainly a threshold effect : to emit x neutrons a certain 
value of the energy must be reached. When the energy reachs the siiitable value, the fission time increases due to 
the stronger competition between fission and evaporation. Since each neutron emitted carries some kinetic energy, 
proportionally to the temperature of the nucleus under consideration (square root of the energy), the period of the 
oscillations increases with excitation energy. However peaks are separated by roughly 20 MeV, which is a too large 
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value compared with the binding energy of a neutron. If we remove the pairing energy from the mass of the nuclei, 
new peaks appear between the peaks shown on figure 3. This succession of peaks with a periode of approximatively 
10 MeV can now be easily understood as a threshold effect. Thus the oscillations and their behavior come from 
the combined effects of energy threshold and pairing energy. We notice also that as we go to higher energy, the 
amplitude of the oscillation decreases. The damping of those oscillations can be explained by the increase of the 
number of steps of evaporation at higher energy. This conduces to a boarding of the energy distribution. Thus the 
threshold effect becomes less sensitive and the amplitude of the oscillation is reduced. At low excitation energy the 
disagreement between the calculations and the experimental data looks more important. This can be attributed to 
the quantum effect that persists at such low temperature of the compound nucleus. This hypothesis is confirmed by 
the experimental observation of asymmetric fission below 50 MeV in this experiment. 




FIG. 3. Fission time of^^^U as a function of the excitation energy. The upper panel represents the calculation done with a 

step of excitation energy of 3 MeV with angular momentum of the compound nucleus J — 1,5 and 9h. The lower panel shows 
an averaging of the calculation on a range of 30 Me V. The curves represent the calculation done with various values of the 
friction factor (3. The symbol represent the experimental data. 

We see that the experimental data present no oscillation, but have very large error bars. The huge size of the error 
bars would be the main argument to explain the absence of oscillations. Actually, if we plot the same calculation as 
in the upper panel but averaged on a range of 30 MeV (range comparable to the mean value of the error bars) shown 
in the lower panel of the figure 4, we can see that the oscillations roughly disappear. 

The agreement between the experimental data and the calculation are reasonably good. The experimental value and 
the calculation are all in agreement within a range of roughly one order of magnitude. The evolution of the calculated 
curves as a function of the friction parameter /? is rather strong when we go from /? = 1 lO^^s"^ to /? = 5 lO^^s"^ 
this is due to the Kramers factor (see equation (37)) as describe in the section 1. It seems that the friction value 
of 5 lO^-'^s"^ is the best to reproduce experimental data. However, we have to be careful on some points to give 
definitive conclusion: information is missing on angular momentum distribution. The reliability of the conclusion 
could probably be improved by having more information on the distribution of angular momentum as a function of 
the excitation energy. 

We should notice that in those calculations the value of the reduce friction coefficient (3 that gives the best agreement 
with the data is 5 lO^^s"^ and in the section IV A the value of /3 is set to 5 10^°s~^. This difference in the suitable 
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value of /3 is a rather delicate problem which goes beyond the scope of this paper. Detailed analysis will be published 
elsewhere. 



V. CONCLUSION 

The KEWPIE code has been designed as a disintegration code dedicated specifically to the study of superheavy 
elements. Its conception is based on the idea of minimising the number of free parameters. So the only free parameter 
that is adjusted in the calcidation is the reduced friction coefficient /3. We have shown that the results of the 
calculation with such a high level of constraint of only one adjustable parameter are rather successful. The quality 
of the agreement is as good as other evaporation codes that are including more adjustable parameters. The results 
of the calculation on the superheavy elements show that the shell correction energy is a key physical quantity for 
understanding and reproducing the residue cross sections. 

The architecture of the program as a full cascade is an innovation which permits us to include the dynamics of the 
reaction. Thus, the evolution during the decay can be calculated. This time dependence permits to calculate fission 
time dynamically. The results of the calculations appear to be in quite a good agreement with the experimental data. 

The KEWPIE program permits to study the disintegration of nucleus and to accommodate both a statistical and 
a dynamical frameworks. Moreover the way of programming allows us to implement new physical effects in the code 
without affecting the whole structure. 
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